﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace Nature.Numerics.Basic
{
    public abstract class Newton
    {
        public abstract double Solve(double x);

        #region Private Fields
        private double m_absTol = 1.0e-5;
        private double m_relTol = 1.0e-5;
        private double m_absRadius = 0.0;
        private double m_relRadius = 1.0e-1;
        #endregion

        #region Properties
        public int MaxNumberOfIterations
        {
            get;
            set;
        }

        public double AbsTol
        {
            get { return m_absTol; }
            set { m_absTol = value; }
        }


        public double RelTol
        {
            get { return m_relTol; }
            set { m_relTol = value; }
        }

        public double AbsRadius
        {
            get { return m_absRadius; }
            set { m_absRadius = value; }
        }

        public double RelRadius
        {
            get { return m_relRadius; }
            set { m_relRadius = value; }
        }

        public int NumberOfIterations
        {
            get;
            protected set;
        }

        #endregion

        public static Newton Create<TFunc>(TFunc function) where TFunc : IFunction
        {
            return new Newton<TFunc>(function);
        }


        public static Newton Create(Func<double, double> f) 
        {
            var genericFunction = new GenericFunction(f);
            return Create<GenericFunction>(genericFunction);
        }

        class GenericFunction : IFunction
        {
            Func<double, double> m_func;
            internal GenericFunction(Func<double, double> f)
            {
                m_func = f;
            }

            #region IFunction Members

            double IFunction.Compute(double variable)
            {
                return m_func(variable);
            }

            #endregion
        }

    
    }


    public class Newton<TFunction> : Newton
        where TFunction : IFunction
    {
        public Newton(TFunction function)
        {
            Function = function;
            MaxNumberOfIterations = 100;
        }

        public TFunction Function
        {
            get;
            private set;
        }

        public override double Solve(double x)
        {
            NumberOfIterations = 0;
            double dx, f, dfdx, del, radius;
            bool solved = false;
            while (false == solved)
            {
                NumberOfIterations++;
                if (NumberOfIterations > MaxNumberOfIterations)
                {
                    throw new System.InvalidOperationException();
                }
                f = Function.Compute(x);
                del = 0.5 * (AbsTol + Math.Abs(RelTol * x));
                dfdx = (Function.Compute(x + del) - f) / del;
                dx = -f / dfdx;
                radius = 0.5 * (AbsRadius + Math.Abs(RelRadius * x));
                double ratio = Math.Abs(dx / radius);
                if (ratio > 1.0)
                {
                    x += dx / ratio;
                    continue;
                }
                x += dx;
                solved =
                    (
                        Math.Abs(dx) < AbsTol ||
                        Math.Abs(x / dx) < RelTol
                    );
            }
            return x;
        }
    }
}
